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Abstract 

We study the propagation of localised disturbances in a turbulent, but momen- 
tarily quiescent and unforced shell model (an approximation of the Navier-Stokes 
equations on a set of exponentially spaced momentum shells). These disturbances 
represent bursts of turbulence travelling down the inertial range, which is thought 
to be responsible for the intermittency observed in turbulence. Starting from the 
GOY shell model, we go to the limit where the distance between succeeding shells 
approaches zero ("the zero spacing limit") and helicity conservation is retained. We 
obtain a discrete field theory which is numerically shown to have pulse solutions 
travelling with constant speed and with unchanged form. We give numerical evi- 
dence that the model might even be exactly integrable, although the continuum 
limit seems to be singular and the pulses show an unusual super exponential decay 
to zero as exp(— const u") when n — > 00, where a is the golden mean. For finite 
momentum shell spacing, we argue that the pulses should accelerate, moving to 
infinity in a finite time. Finally we show that the maximal Lyapunov exponent of 
the GOY model approaches zero in this limit. 
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Shell models of turbulence have increasingly been used as a laboratory for 
testing hypotheses about the statistical nature of Navier-Stokes turbulence. 
Even though they are derived by heuristic arguments they are surprisingly 
adept at reproducing the statistical properties of high Reynolds number tur- 
bulence. The reasons behind the apparent success of the shell models still 
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remains unclear. In the present paper an attempt will be made to examine 
the detailed behaviour of one particular shell model - the so-called GOY model 
by Gledzer-Ohkitani-Yamada [HO] in the limit where the ration between the 
wave number of successive shells approach unity - which we shall refer to as 
the zero-spacing limit. 

One of the most important problems in turbulence is to understand inter- 
mittency - the very lumpy and irregular intensity of turbulent motion. Even 
in strongly turbulent states each region experiences periods of almost com- 
plete quiescence interrupted by rapid burst events, and this is believed to be 
the origin of the corrections to Kolmogorov scaling seen experimentally, and 
predicted rather accurately by the shell models [3]. Single burst motion su- 
perimposed on the Kolmogorov spectrum (which emerges as a fixed point of 
the shell model) have been studied in shell models by a number of authors 
[11-14]. In this paper we also consider single burst motion but, in contrast 
to the previous work, we do not superimpose a Kolmogorov background. It 
was shown recently, that prior to a burst in the GOY shell model, the veloc- 
ity field tends toward the "trivial" quiescent fixed point [15]. In this paper we 
thus study localised disturbances propagating through the quiescent, unforced 
system. One should keep in mind in the following, that the shell models are 
defined in /c-space and that the pulses we describe are likewise in /c-space. 
Thus a pulse travelling to smaller k should correspond roughly to a spatial 
region of growing size, but we shall not go into the detailed conversion of our 
results to real space, since this is not even a completely well defined process 
in the context of (1-d) shell models. 

The paper is structured as follows: First a short introduction to the GOY 
model is given in Section 1. It was conjectured by Parisi in a yet unpublished 
paper [12] that the correlation functions of the GOY model (also called the 
structure functions) might be calculated from the existence of soliton-like so- 
lutions on top of the Kolmogorov fixed point u{k) oc k~^^^. For the study of 
the pulses he proposed a continuum limit of the GOY model (the "Parisi equa- 
tion"). Around the quiescent fixed point this equation does not form pulses, 
but rather forms shocks. The Parisi equation around the quiescent fixed point 
is studied in [16] and the results are briefly reviewed in section 2. Pulses in 
the quiescent state are, however, seen in the zero-spacing limit of the GOY 
model. The equation governing the zero-spacing limit is derived in section 3 
and the qualitative behaviour of our model is illustrated through numerical 
simulations. In section 4 a continuum version of the zero-spacing limit is pro- 
posed, and a solution for a pulse is derived. The asymptotic behaviour of a 
pulse is is treated in the case where the background field is zero (section 5) 
and non-zero (section 6). In section 7 the acceleration of a pulse which is seen 
away from the zero-spacing limit is examined. Finally the Lyapunov exponent 
of the GOY model is calculated as a function of the shell spacing, and it is 
shown that the model is not chaotic in the limit of the shell spacing going 
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towards zero (Section 8) 



1 The GOY model 



Shell models are formed by a truncation approximation (in A;-space) of the 
Navier-Stokes equations [10]. The most well-studied model is the "GOY" 
model of Gledzer-Ohkitani-Yamada. 

This model yields corrections to the Kolmogorov theory [3] in good agree- 
ment with experiments [17,18]. For the "GOY" shell model, wave-number 
space is divided into separated shells each characterised by a wave-number 
kn = r" ko with n = 1, - ■ ■ ,N. Each shell is assigned a complex amplitude 
describing the typical velocity difference over a scale = By assuming 

interactions only among nearest and next nearest neighbour shells and phase 
space volume conservation one arrives at the following evolution equations [2] 



(^ + O - 

h c 

« ^n(an 1in+l1*n+2 + — ^^n-l^^n+l + + (1) 

with boundary conditions = h^q = C\ = C2 = a^-i = ajv = 0. / is an 
external, constant forcing, and u is the viscosity. 

An important property of the GOY model is that it is possible to incorporate 
the conservations laws of turbulence as found in the Navier-Stokes equations. 
The first conservation law that needs to be satisfied is the conservation of 
energy. The energy in the GOY model is given by 

E = El^nl' (2) 
n 



which should be conserved in the limits of no forcing and vanishing viscosity. 
This leads to the following relation between the coefficients 

On + K+1 + Cn+2 = (3) 



The constraint still leaves a free parameter 6 so that one can set a„ = 1, bn+i = 
—S, Cn+2 = ~(1 ~ ^) [9] (the value of is fixed by the time scale). Kadanoff 
et al. [6] observed that also a "helicity" invariant exists on the form 

H = ^{-irkn\Un\' (4) 
n 
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Conservation of this quantity leads to an additional constraint on the param- 
eter bn 

bn = ^-l or 6 = 1-^. (5) 

In the most studied case where r = 2 one obtains the canonical values 6„ = 
c„ = |. In this paper we are going to focus on the limiting case r — > 1 and 
according to Eq. (5) 6„ — > in order to satisfy helicity conservation. The set 
(1) oi N coupled ordinary differential equations can be numerically integrated 
by standard techniques. 



2 Some results on the Parisi continuum limit 



In a previous paper [16] we studied a number of features of the Parisi [12] 
continuum model. In particular, a "burst" event (the evolution of an initial 
condition having support only in a finite interval of /c-space) was studied. The 
initial burst splits up into a right and a left moving part and on the front 
of each part a shock is formed. It was shown that the continuum equation 
without forcing or dissipation can be explicitly written in characteristic form 
and that the right and left moving parts can be solved exactly. When this 
is supplemented by the appropriate shock condition it is possible to find the 
asymptotic form of the burst. 

In this section we shall give a brief review of some of the main features pre- 
viously obtained. The continuum Parisi equation is obtained from the GOY 
model by letting the shell spacing approach zero, and expanding to the first 
order. By writing the distance between the shells as r = 1 -|- e with e <^ 1, 
we have kn ~ exp(ne), so with n ~ const. /e a continuous range of values is 
obtained for the variable k. 

To proceed we use a Taylor expansion of the type 

du. 



Un+i(t) ^u{kn+i,t) = u{\nkn+i,t) ft! 'u(nln(l + e),t) + ln(l + e) 
^u{k,t) + ek 



n 



dink 

du{k, t) 



dk 

— u(k, t) + ek ^^[^^ , with u(k, t) = u(kn, t) — Un(t), (6) 
ok 

and similarly for Un+2,Un-i, and Un-2- To first order one obtains [12]: 

u; + iyk\* = -ik{2 - S) (u^ + Skuuk) , (7) 
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where Uk = du/dk etc.. By rescaling time with 2 — 5 we get the Parisi equation: 

ul + Sikuuk ~ —iku^, (8) 

for the case of no forcing or viscosity. Here subscripts denote differentiation. 
Sphtting into real and imaginary parts, u — a + ib, we get 

at — Skabk — Skabk = 2kab (9) 

and 

bt-3kaak + 3kbbk^ k{a^ -b'^). (10) 

In numerical simulations we observed that a short while after the initialisation, 
shocks are formed in both directions. As was seen from the equations, the real 
part is only able to move to the left, while the imaginary part moves both ways. 
Setting a = in the original equations, we get just a single real equation for 
the inertial range: 

bt + 3k%bk^ -kb"^ (11) 



For the left-moving part of the pulse, the real and the imaginary parts are 
seen to become proportional. Inserting b = Ca into the Parisi equation gives 
us two possibilities: a = (which we already treated) and a = ±\/36. For 
a — y/Sb the equation for the imaginary part becomes: 

bt - 6k%bk = 2kb'^. (12) 



It is possible to explain the observed splitting of the pulse by writing (8) in 
characteristic form. We shall not give the details, since they can be found in the 
previous paper [16]. The starting point is to rewrite (8) in a form resembling 
the the Burgers equation 

— ivVx = 0, (13) 



where u = k^^/'^v and k = (2x)^"^/^. Of course, there are three important 
differences to the Burgers equation: Here we have an equation in fc— space, v 
is complex and the conserved quantity is the energy E — J{2x)^^'^v'^. Writing 
V — re^^l^ the Riemann invariants J± can be found, and the solution can be 
expressed in terms of these 



(14) 
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and 



sin^ = 



J|- J! 



(15) 



The formulation in terms of characteristics can now be used to understand the 
sphtting of the pulse found in Section 3.1. The Riemann invariants J± have 
the property that they are constant along the curve: 



The important point is that, since r is nonnegative, one family of characteris- 
tics ( J+) can never move to the left while the other one (J_) can never move 
to the right. If the initial condition has support in a limited region of x, say 
the same is true of J±. They both vanish (in the initial condition) 
outside of this interval. For times t > we compute the field values by finding 
where, on the x-axis (i.e., t = 0), the characteristics going through the point 
{x,t) emanate. Now, if a: > x+, the J_-characteristic going through this point 
must emanate from some Xq > x+ and thus J^{x,t) = J_(xo,0) = 0. This 
means that either r = (which makes the entire field vanish) or sin^ = 1 
which means that 6 has the constant value 7r/2, which implies that a = 0, in 
agreement with the numerical simulations. 

The splitting of the pulse into a right-moving part with a = and a left- 
moving part with a = \/^b makes it possible to give a complete solution for 
a single burst event. For the details, we refer to [16]. Each pulse develops a 
shock which can be followed by imposing conservation of energy across the 
shock. 

The asymptotical solutions for the two pulses are: 



both decaying in time. Thus a burst created on top of the zero solution not 
only travels down the inertial range until it is dissipated by viscosity (here, 
by the continuum), it also travels upwards, to the smallest k. The burst does 
not remain localised in /c-space, but is distributed over the /c-range, and then 
decays. It is interesting to note that the initial disturbance creates both a 
forward and an inverse cascade. 



dx 
'dt 



±r{x,t). 



(16) 




(17) 
(18) 
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Fig. 1. The development of a initial small bump in the GOY model with r = 1.1 
and r = 1.0. Shown is |n„|. 

3 The zero-spacing limit 

When we study the GOY model for r > 1 we observe important deviations 
from the behaviour of the Parisi equation described in the previous section. 
We focus on the motion on a single pulse in the case where both the forcing 
and the viscosity is zero, i.e., u = f = and the two cases r = 1.0 and r = 1.1 
(Fig. 1). the initial condition is a localised bump: 

Un — Uo{l + i) {cos[(n — no)Ko\ + 1} for x e [no^o — ^o'^o + '^], (19) 

where uq is the center of the bump, 27:k,q is the length of the bump and ^uq 
is the height of the bump. 

The initial bump splits up into a left and a right moving part with the same 
phases as in the Parisi equation. However, a shock is not formed, rather the 
initial bump splits into a number of individual pulses. In the case r — 1.10, the 
pulse to the right accelerates whereas the pulse moving to the left decelerates, 
similar to the shocks in the Parisi equation. In the case r = 1, this does not 
happen and both the right- and left-moving pulses propagate with a constant 
velocity. In both cases, however, the shapes stays invariant as time progresses, 
completely similar to the behaviour of a soliton in an integrable system, like 
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the KdV equation. 



To be able to observe the pulses in more detail, we perform the limit r — > 1 in 
a different way as is done when the Parisi equation is derived. We shall show 
that this leads to an interesting model containing pulses moving with constant 
velocity and which might even be a new example of a completely integrable 
discrete field theory. 

The starting point is the GOY model 

^~dt ^ """"^ ^ -lkn{Un+lUn+1 " -Un-lUn+1 ^Un-lUn-2) (20) 

where A;„ = fcor". We now introduce the variable w„ = kn-iUn, and in terms 
of this variable we get 

+ ^klWn)* = -i(r~^ty„+iWn+2 - 5Wn-lWn+l " (1 " 5)r'^Wn-lWn-2){^l) 

Using this substitution we get rid of the explicit dependence on A;„ on the right 
hand side of Eq. (21). 

Now we take the limit r — > 1 and at the same time wc invoke helicity conser- 
vation, which means that we must choose 5=1 — 1/r^O. Finally, we are 
interested in very large Reynolds numbers, so we let u ^ and then we arrive 
at the model 



The middle term in the GOY model has disappeared and we are left with the 
two terms directly responsible for the propagation of the pulses to the left and 
the right. 

This model has turned out to have very interesting properties. An initial bump 
on zero background {w = 0) splits into a scries of left and right moving pulses 
with constant velocity as already shown in Fig. 1. Similar to the case for the 
Parisi equation (8) the phases of the pulses only have certain fixed values. Let 



dt 



:) 



(22) 



(23) 



then 



dR 



■n 




Rn-lRn-2)e'^'^-''/'^ 



(24) 



dt 



■n+2 — 
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Fig. 2. Time-space plots of the imaginary part of the transformed GOY model Eq. 
(26), with a cosine-bump as initial condition (19). The boundary conditions are 
made periodic, to make it possible to see the interaction of the pulses. The right 
picture is a zoom on an interaction between two pulses, where the time-delay in the 
collision process is clearly visible 

and to have a solution with constant phase, the exponential factor has to be 
real. This means that sin(30 — 7r/2) = or 







2p+ 1 



-tt; 



(25) 
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precisely the same values as for the Parisi equation [16]. For these values of 
(f), the exponential factor in (24) becomes cos(30 — vr/2) = cos{pn) = (— 1)", 
giving rise to the left (+1) and right (-1) moving fields respectively. This 
behaviour is in fact an intrinsic property of the inviscid GOY model, and 
stem from the quadratic terms and the complex conjugation. This property 
of selecting specific phases might explain the period-three organisation of the 
pulses in the forced GOY model as noted by Okkels & Jensen [15]. 

The simplest case, on which we shall concentrate below, is p = 1; a completely 
imaginary field moving to the right. We thus let Wn — iwn and dropping the 
tildes we get the real equation 

-T- = -Wn+lWn+2 + W„_iW„_2. (26) 



The behaviour of (26) using a cosine-bump as initial condition (19) is shown 
in Fig. 2. The initial bump splits into a number of pulses with the same shape 
but different heights. The velocity of each pulse was measured to be w = 2h 
where h is the height of the pulse. Two pulses can pass through each other 
after a time delay produced by the collision (Fig. 2, right). These facts seem 
to indicate that the model (26) is exactly integrable - although, as we shall 
see later, probably not in terms of analytic functions. 

One has to keep in mind that the system is discrete and therefore the shape 
of a pulse deforms slightly in time. Despite this, one can actually obtain a 
very precise shape function for the pulse, simply by measuring the field w 
continuously in time at a given point in the lattice (Fig. 3) 

The most characteristic feature of Eq. (22) compared to the original GOY 
model, is the disappearance of the middle term. It should be noted that other 
variants of shell models display another limit. An example is the hclicity- 
conserving shell models [19] where in the r = 1 limit the right term disappears, 
while the middle term is preserved. This makes the equation qualitatively 
different. 



4 Pulses in the continuum limit 



We shall now show that the continuum hmit of (26) is non-trivial and thus that 
the form of the pulses might be unusual. By expanding w in a Taylor-series 
and retaining, on the right-hand-side of (26) terms w^^w^'^^ with g < 3 we 
get the continuum limit 

-wt^ 6wwx + QwxWxx + ^wwxxx (27) 
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Fig. 3. To the left is shown the shape of the pulse. The circles is the pulse extracted 
directly (at t = 210 in Fig. 2) and the line is the pulse extracted using a time series 
from one point. On the figure to the right the horisontal axis is z = o"", where 
a = 1.618... is the golden mean. Thus the tail of the pulse fits well to the form 



where again the subscripts denote differentiation. 

Using the identity {w^)xxx — ^w^w^x + '^wWxxx the continuum model can be 
written: 

-Wt = 6wWx + {w'^)xxx - WWxxx (28) 



To study the asymptotics, we shall start by omitting the last term. The re- 
sulting equation 

-wt^ 6wwx + {w'^)xxx (29) 



strongly resembles the Korteweg de Vries (KdV) equation. We shall refer to 
it as the Quadratic Korteweg de Vries (QKdV) equation. 
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4-1 Solitons in the KdV equation 



Let us briefly review how the sohtons of the KdV equations are found. We 
look for pulse solutions of the KdV equation 

Wt = WWx + Wxxx (30) 



which are only a function oi x ^ x -\- vt, so that (30) becomes 

vw' = ww' + w'" (31) 



Now we can integrate once to get 



w" -VW + ^w^ = Ci (32) 



Since we are interested in a pulse where the LHS approaches as x — > ±oo 
we must take C\ — 0. The resulting equation is that of a Newtonian particle 
moving in time x in the potential 

[ \ 11 

Y{w) = / {-vz + -z^)dz = --vw^ + -w^ (33) 

J Zi Zi \j 



Energy conservation has the form 

\{w'f^V{w)=C^ (34) 



Again the boundary conditions dictate that C2 = and we can solve (34) as 

dz 

I (35) 
.0 V-2^(^) 



This leads to the sech'^ for w{x) but for our purposes the most important 
point is that the soliton corresponds to the homoclinic orbit starting from 

the potential hill-top w = Oata;^— 00 and returning back at x —00. 
The quadratic nature of the hill-top means that the integral for x{w) diverges 
logarithmically for w — > and thus the soliton extends to infinity with an 
exponential tail in both directions. 
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4-2 Solitons in the QKdV equation on a vanishing background 

We now look for a pulse solutions of (29) and let x — > x — vt, so that (29) 
becomes 

vw' = 6ww' + {w^Y" (36) 

Again we can integrate once to get 

{wy -VW + Sw'' = Ci (37) 

and, if at some (finite or infinite) point the left hand side vanishes, we must 
take Ci = 0. This means that we are looking for pulses with w — >^ c = 0. To 
make a Newtonian problem out of this we now introduce the variable y — w^. 
Then we get 

y"-v^ + 3y^C, (38) 



which is potential motion in the potential 



V{y) = -Iby'^' + \v' (39) 



where b — Av/9 and the total energy (34) must again be 0. Thus the solution 

can again be written as (35), but this time the integral converges at ty = so 
we can write explicitly for the solution which vanishes at a; = 0: 



x(w] 



W 

ds 2 f dz 2 . b - 2w , 

- / , = COS"^ ; (40) 

J s^bz - a/3 b ^ ^ 



-2V{s) VsJ Vbz - z^ 



which can be inverted as 



b VSx 2v ^/3x 

«; = -(l-cos^) = y(l-cos^) (41) 



We verify that w(0) = but now the solution has support only on the interval 
[0, 47r]. Outside of that it is identically zero which means that there is a shock 
in Wj::^ at those points. The strength of this shock will presumably diminish 
when higher order approximations to (26) is used. But one can easily see that 
the term ww^xx which we omitted from (28) will not change this fact (although 
it will probably alter the simple solution (41)). Closely to the right of a; = 0, 
w{x) ~ ^x"^. Thus Wx and {w'^)xxx both go linearly in x (and balance), whereas 
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{w'^)x and wWxxx both go like and do not alter the asymptotics for x — > 0"^. 
This will be verified explicitly in the next section. 



4-3 Higher approximations. 



It is possible to take the term wWxxx, which we omitted in (27) explicitly into 
account [20] . Multiplying (27) with 2w on both sides, we can rewrite this as: 

- {w% = + 6{w''wxx)x (42) 



If we again assume an uniform propagation with constant velocity v, then 
w{x, t) is only a function of a; — vt, and the equation above reduces to 

v{w''y = A{w^y + 6{w^w"y (43) 



which can immediately be integrated to yield: 

w^{-v + Aw + 6w") = constant (44) 

The non-trivial solution must therefore satisfy 

6w" + 4w = v (45) 

with solution in the form of the "compact" soliton 

V [2 
w ^ ^{l - cos ^ -x) (46) 

analogously to (41) but no solutions decaying to at infinity 

This approach is tied to the existence of a conserved energy given by: E = 
J dxw{x, ty, which is conserved only if 2wwt can be written as a total deriva- 
tive in x: 

2wwt^ F{w,Wx,..)x (47) 



But if we furthermore assume a dependence only on x — vt, then the above 
equation reduces to an equation of total derivatives: 

-v{w^y = F{w,w',..y (48) 
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which can be integrated to give: 

F{w, w', ..) + vw'^ = constant. 



(49) 



To have finite-energy sohtons, w must go to zero at infinity. If furthermore 
F{w, w', ..) goes to zero as w goes to zero, the constant must be zero. We can 
then divide by and the equation reduces to 

Ei^^i^iA ^ (50) 

Keeping terms up to 5'th order in the derivatives in equation (26) in the 
continuum hmit, we get: 

11 3 

xxxx ~l~ '^^xx^xxx- 

(51) 

Miraculously, this can be rewritten as 

/ 11 2 2 4 \ 

- {w\ = \^{W^) + ^{w'^Wxx) + l^w'^Wxxxx + -^WWxWxxx " -^wlw^x + (^2) 



so that 



11 2 2 4 

F(w, w\ ..) = -4w^ - Qw'^Wxx - —w'^Wxxxx - -WWxWxxx + T'^lWxx " -ww^^(53) 

zU o 00 



Equation (50) is then: 



v = Aw + Qw" H w"" + 



11 , 2w'w"' + 2w"'^ 2w'^w" 



20 5 w 5 

^ „ 11 „„ 2 /w'V'V 1 
20 b \ w J w' 

which might have non-trivial soliton-solutions as long as the last term equals 
V at large x — vt. So for large x — vt, we must have: 

-vw — 55 

2 \ w ^ ^ 



or, by integration: 



5 w'^w" 

-vw — C 56 

2 w ^ ^ 
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where C is some integration-constant. Now multiply this with ww': 

^v{w^y + 2C{w^y ^{w"")' (57) 



and integrate again 



^vw^ + 2Cw^ = w'^ (58) 



where no extra constant arises, because of the boundary-condition at infinity. 
Thus 

/ID?; \^/^ 

yj'= (^yj5 + 2CwA (59) 



If C 7^ 0, the 2nd term dominates asymptotically, and we get w{x, t) ~ {x—vtY 
for large x — vt. On the other hand if C = then w{x, t) ^ (x — vt)^ for large 
X — vt. In both cases w{x, t) does not decay to zero as a; — goes to infinity. 
So again there are no soliton-solutions, that decay to zero at infinity. 



4-4 Solitons in the QKdV with c > 0. 



The strange appearance of the solitons is closely linked to the fact that w 
decays to zero asymptotically. If w — > c > for x — > ±oo, the pulse will decay 
exponentially to this value far away just as in the KdV case. The difference is 
that the constant Ci now no longer vanishes. Instead Ci = —vc -\- (? and the 
potential acquires the form 

y{v) = -\^"' + \y' - C,y (60) 



Now, for i> > and small c > 0, Ci —vc < so the potential acquires a 
small hill top aX w = c. This hill top is non singular (albeit in a very small 
region) and we get back to the usual situation (like KdV) where the soliton 
spreads over the entire region with an exponential tail. Near w = c we get 
from (35) 

y 

f '^■'^ 1 y — c 

xiy) — xiyo) = / , = , =log (61) 
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where V"{c) — — ^. Thus the solution decays as 

y{x) = w\x)=c+{yo-c)e'''' (62) 

where a — —V"{c) ~ (y^)~^. Thus the decay length ^ ~ ^/c. 
5 Asymptotics of solitons in the discrete model 



The results obtained in the continuum limit in the previous section are not 
born out by simulations. The pulses do not seem to have compact support, 
but decay far away. The decay in not exponential, but seems to more rapid. 
We shall now give an estimate of this decay and subsequently try to refine it. 

Let us approximate the derivative in (26) as w'{x) ^ w{x + l) — w{x). Further, 
expecting a rapid decay, we shall neglect w{x + 1) and w{x + 2) compared to 
w{x), w{x — 1) and w{x — 2). Thus we simplify (26) to 

w{x) — w{x — l)w{x — 2) (63) 



Now let ({x) — log(x) whereby we find 

C(x) = C(x-l) + C(a;-2) (64) 

This recursion relation is identical to the one for the Fibonacci numbers, and 
in general ({x) will grow like a^. To determine a, assume that 

C{x) = acr^ (65) 

which should be vahd at large x. Inserting into (64) gives us 

- a - 1 = (66) 



with the (positive) solution 



a = i(l + V5) (67) 



the so-called golden mean. This gives the extremely rapid decay 

w{x) « e-"'"" (68) 
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with some positive constant a. This rather unusual behaviour is extremely 
well represented by the numerical solution (Fig. 3). 

To refine this, assume that 

w{x) = /(a;)e«^' (69) 

Inserting this into (26) and canceUing a factor of e"^"" (using (67)) we get 

af{x)a^ log a - fix) = f{x - l)f{x - 2) - f{x + l)f{x + 2)6--^-^ (70) 

where d = a + a"^ — a~^ ^ 3.2 > 0. If the function f{x) is "reasonable" 

the last term can be dropped and we get 

af{x)a^ loga - f{x) = f{x - l)f{x - 2) (71) 

and we can get the asymptotic behaviour / ~ cx^ by matching the first term on 
the left hand side to the right hand side. Thus f{x) ca^ with c = aa^ log a 
and the solution takes the form 

w{x) = aa(^+2) logae""' = alogae("+2)^°s'^+"^' (72) 

Note that a is still not determined and we presumably need to solve for the 
whole soliton structure to get this constant. 



6 Expansion for non-zero background 

As seen in section 5 pulses converging to cq > decay exponentially. This is 
born out by out numerics and we now give a method to obtain the asymptotic 
behaviour of a pulse solution. 

We consider 

Wn = Wn+lWn+2 - Wn-lWn-2- (73) 

For the field Wn{t) we make the expansion 

s=oo 

'^n{t) = X] Csex.-p{-s{kn - vt)). (74) 

s=0 
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0.5 1 1.5 



Fig. 4. An example of an "exact" solution of a pulse using the expansion Eq. (74) 
with the coefhcicnts obtained by the recursion relation Eq. (76) with k = 2, cq = 0.5 
and Cl = 2. Note that Wn asymptotically approaches cq. 

Hence we have 



oo oo I 

V Ici exp{—l{kn — vt)) = — 2 y^exp(— /(/cn — vt)) y^CsQ_g sinh(2ZA: — sk). (75) 



1 



1=1 



s=0 



Thus we have the recursion relation 



ci{vl - 2co sinh(2ZA;) - 2co sinh(ZA;)) = -2 X] CsCi-s sinh(2ZA; - sk). (76) 



s=l 



Although the recursion relation is non-linear, it turns out that it can be solved 
successively in terms of the two constants Cq and Ci. For Z = 1 we get 



V = 2co(sinh k + sinh 2k) , 



(77) 



giving a relation between the constant background cq, the width k, and the 
velocity v. For Z = 2 we obtain 1 



C2 = cl 



2 sinh 3A; 



2v — 2co sinh 2k — 2cq sinh Ak 



(78) 
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Similarly, for Z = 3 we get 

^ 3 (sinh 5k + sinh 4A;) sinh 3A; 

'^^ (3^; — 2co sinh 6k — 2co sinh 3k) {2v — 2co sinh Ak — 2cq sinh 2k) 



It is straightforward, but tedious, to continue this calculation (we have com- 
puted C4 and C5). 

To show that (74) is a solution, we need to show that the coefficients q ob- 
tained from the recursion relation (76) are such that the sum in (74) is either 
convergent or Borel summable. We have not succeeded in this in general. How- 
ever, it is possible to show convergence when the width k is large, and also 
presumably when k is small. Let us consider k to be large. The velocity then 
becomes 

V^coe^'' + O{e''). (80) 

^Prom the recursion relation (76) we get a solution of the form 

ci ^ cle-'^^'-^Vco^'ll + ^(e-'')). (81) 



The corresponding asymptotic solution can be found by performing the sum 
in (74), 

~ 1 _ (ci/co)e-fc('^+i)+'=o*exp2fc- (S2) 

Numerically, one can easily generate an "exact" solution by Eq. (74) using the 
recursion relation Eq. (76). Of course the solution exists only if the coefficients 
Cs converge. We ohave observed that the radius of convergence in parameter 
space has a "hole" around the point Cq = 0, Ci = 0. This is not surprising 
because this expansion is only valid on non-zero background and from the 
structure of the recursion relations it is clear that Cq = Ci — leads to a 
trivial case where = for all s. Away, from the open set around co — ci — 
we find a very fast convergence of the series. An example of the obtained 
solution is shown in figure 4. 

Note that formally we can apply this expansion also to the case treated in 
section 6, where w — > for r — > 00: At first sight this seems to require Cq — 0, 
which doesn't work. However one has to keep in mind that the results of the 
previous section shows that in this case k must be negative, indeed k = —loga 
and this means that we have an expansion in terms of growing exponentials. 
The background value c is therefore not given simply by cq in this limit, but 
all higher terms will contribute to it. 
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7 Acceleration of solitons for r > 1 



When r is slightly above unity, the main difference between our zero-spacing 
model and the real GOY model seems to be the acceleration of the bursts 
seen e.g. in Fig. 1. Letting r = 1 + e in (21) and retaining helicity conservation 
(5 = 1 — 1/r ~ e) we see that our inviscid model (22), in lowest order of e, 
becomes 

dw* 

= -i{wn+iWn+2 " w„_iW„_2 - 4e«;^) (83) 

For purely imaginary w this is 
dwn 



dt 



-Wn+iWn+2 + w„_iW„_2 + 4e«;„ (84) 



In the continuum limit, neglecting all higher order derivatives we get 

— = -&ww^ + Aew^ (85) 
at 



from which we clearly see the acceleration caused by the last term. Note that 
this equation preserves energy conservation in the form E — J u^dx — const, 
but now we have to remember the transformation leading to (21) , i.e. u — wk, 
where k — ko{l + e)" ~ /cqG^^. Thus (85) can be written 

^{w^e-^'^) = -^{4w^e-^'^) (86) 



which shows that the appropriate shock-condition is 



where V is the velocity of the shock and [f{x)\ denotes the discontinuity of 
the quantity x accross the shock. Since w vanish on one side of the shock we 
get 

V ^Aw (88) 



The characteristic equations for (85) are 



dx „ 

— = 6«; (89) 
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and 



(90) 



The latter can be integrated as 



w = (91) 
1 - 4:ewot ^ ' 



where is the initial field taken as a function of the initial coordinate xq. 
The second characteristic equation can now be integrated as 

3 

X = xo - — log(l - 46^0^) (92) 



Thus each characteristic ends in a singularity at time t — l/4ewo(a;o). Never 
the less it is perhaps not obvious that the shocks also have to end in a singular- 
ity at a finite time. To investigate this, we must investigate the shock condition 
together with the solutions (91) and (92). The shock velocity is y = dXf,/dt, 
where is the edge of the shock. The natural variable is however xq and we 
thus express the shock-condition as 

dx dx dxn dx , „, 

dt dxo dt dt ^ ^ 



On performing the difi^crentiations of (92) with xq and t and inserting the 
solution (91) for w we get 

+ <«^* 



To make further progress we must know the form of Wq. In our simulations 
of bursts, we took distributions like (41), which are strictly zero outside the 
interval [a;i,a;2] and approach the ends of the interval with zero slope. Thus, 
close to xi (which must dominate at large times) the form will be 

Wq a{xo — Xi)"^ (95) 



With this assumption (94) becomes: 



(l + 4ay(3-ey)i)^ = -2ay2 (96) 
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where y = xo ~ xi. When t ^ oo and ?/ — >■ the parenthesis on the left 
hand side must be dominated either by the term proportional to t or by the 
constant. If we assume that 4a7/(3 — ey)t — > we get the simple equation 



with the solution 

y^{c + 2at)~^ (98) 

where c = l/yo- With this solution, however, 4aj/(3 — ej/)t — > 6 in contradiction 
to the assumption. Therefore we must assume that yt does not decay to zero 
and the dominant terms are now 

3t| = -y (99) 

which gives the decay y ~ t— 1/3. This solution doesn't work either, since 
the solution (91) only makes sense asymptotically if ty^ — * as t — oo. We 
therefore conclude that, within the approximation leading to the continuum 
limit (85), the position of the shock will diverge to infinity at a finite time. It 
is not clear whether this behaviour is seen in the discrete equation as well(84). 

It is then clear that a pulse reach the dissipative range in finite time. This 
behaviour is qualitative the same as was found for the Parisi equation. 



8 Lyapunov exponents 



In this section we study the Lyapunov exponents of the GOY model in the 
limit of r — > 1. Since we aim at keeping the size of the inertial range fixed, 
the number of shells must also be varied as a function of r, if the "external" 
parameters, i.e., viscosity and forcing are kept fixed at z/ = 10~^, ko = 2~*^, / = 
(l + i)0.005. For r = 2 these parameters correspond to = 19 shells; therefore 
we use a dependence between and r as: 

Ar = 8 + [lli^] (100) 
log(r) 



where the brackets stand for integer value. We wish to quantify the "strength" 
of the intermittency as a function of the effective shell spacing r. One way to 
do this is by means of the maximal Lyapunov exponent [21]. To compute the 
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maximal Lyapunov exponent in the GOY model, we introduce the notation 
U = {Re{ui), Im{ui), ■ ■ ■ , Re{uiy), Im{uN)) and Fj = dUi/dt and consider 
the linear variational equations 

dz- 

=EA,-^, ^ = l,...,2N (101) 



for the time evolution of an infinitesimal increment z = 5U, where 

A^, = dFr^/dU^ (102) 

is the Jacobian matrix of Eqs. (1). The solution for the tangent vector z can 
thus be formally written as z(t2) = M(ti, t2) ■z(ti), with M = exp J^*^ A(r)(iT. 
A generic tangent vector z(t) is projected by the evolution along the eigenvec- 
tor e^^) , belonging to the maximum Lyapunov exponent, i.e. z(t) = |z(0) | e*^^) exp(A 
leading to 

Ai = lim -Inj^ (103) 
t-*oot |z(0)| ^ ' 

where z(0) is the initial tangent vector. 

Practically, Eqs. (1, 101) are integrated simultaneously over a certain time 
5t, starting with a normalised tangent vector in a random direction, z(0). 
The increment over time 5t in the length of the tangent vector is then 5zi — 
|z(5i)|/|z(0)|. Next, the tangent vector is normalised z{5t) — z{5t) /\z{5t)\ 
and this vector is used as a seed for a new integration over the time 5t i.e. 
propagated forward to t = 25t. Generalising this argument we obtain the i'th 
increment 5zi = \z{idt)\/\z{{i — l)6t) \ and the maximal Lyapunov exponent is 
given by (where we now set A = Ai): 

A= lim — 104 

1=1 

We have followed two paths for the values of the parameters an. bn, Cn- In 
the first, the constraint Eq. (5) is applied, and in the other we keep the the 
parameters fixed at the canonical values a„ = 1, 6„ = c„ = — |. The results of 
the numerical simulations are shown in Fig. 5. In both cases we observe that 
the maximal Lyapunov exponent appears to vanish in the limit r — > 1 (the 
minimal value of r in the plot is 1.004, corresponding to ~ 1000 shells). The 
lower curve is the one related to the conservation of both energy and helicity. 
The maximal Lyapunov exponent remains in this case nearly constant in a 
large interval of r- values and then finally drops to zero. In the case where only 
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Fig. 5. The maximum Lyapunov exponent \,nax for a case with conservation of 
hehcity (i.e. 5 = 1 — 1/r) and where 6 is fixed at 1/2. The sohd hne is an interpolation 
to guide the eye. 

energy is conserved (the upper curves) there is a maximum around r = 1.1 
after which also this decreases towards zero. Based on these numerical results 
we conclude that the GOY model does not exhibit chaotic dynamics in the 
zero-spacing limit. Similar behavior has been observed recently in a shell model 
containing two fields [22] . Here it has been observed [23] that the intermittency 
corrections disappear in the limit r — > 1. 



9 Conclusions 



The main objective of this paper is to study the dynamical behavior of single 
burst motion on a vanishing background in various shell models. The motion 
of single bursts is important for the understanding of intermittent behavior 
in turbulence. Indeed it is believed that the presence of intermittency causes 
the corrections to the classical Kolmogorov theory which manifest itself as 
multiscaling of higher order structure functions. In the present paper we have 
studied the motion of a single burst in three different versions of shell models: 

(1) The standard GOY model when the shell spacing approches, but is dif- 
ferent, from zero; 

(2) The Parisi continuum limit of the GOY model; 

(3) The zero-spacing of the GOY model where helicity conservation is kept. 
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In all cases we observe that an initial disturbance splits up in a left- and a 
right-moving part. In case 1) the pulse retains its shape and whereas the right 
moving part accelerates, the left-moving part decellerates. In case 2) an initial 
disturbance splits into left- and right-moving shocks. The right-moving shock 
accelerates and its position diverges to infinity in a finite time. Finally, in case 
3) the shape retains its "solitary" form and does not turn into a shock. Here 
the velocity also remains constant and the solitons pass through each other 
with a slight phase shift but without changing their shapes, completely as 
found in integrable field theories. 

We have analytically calculated the motion of the shocks in case 2) by means 
of characteristics. In case 3) we have analytically estimated the asymptotic 
shape of the solitary wave and find that it decreases super exponentially with 
a decay constant given by the golden mean. On a non-vanishing background 
we have on the other hand found that the decay of the soliton is a stan- 
dard exponential. From numerical simulations in case 1), we have estimated 
the maximal Lyapunov exponent and found that it vanishes in the limit of 
zero-spacing, although the the motion remains chaotic for finite spacing. This 
results is expected since, as we have shown earlier, the continuum Parisi limit 
resembles the Burgers equation and gas theory and is in fact exacly soluble 
by the method of characteristics. 

We hope that our findings could be a starting point for a more detailed un- 
derstanding of structure of intermittency, in particular of how the presence of 
isolated bursts give rise to corrections to classical scaling. 
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